{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.7.3"
    },
    "colab": {
      "name": "example.ipynb",
      "provenance": [],
      "collapsed_sections": []
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "colab_type": "text",
        "id": "ywj1qEt3-lLC"
      },
      "source": [
        "##### Copyright 2020 Google LLC.\n",
        "\n",
        "\n",
        "Licensed under the Apache License, Version 2.0 (the \"License\");\n",
        "you may not use this file except in compliance with the License.\n",
        "You may obtain a copy of the License at\n",
        "\n",
        "    https://www.apache.org/licenses/LICENSE-2.0\n",
        "\n",
        "Unless required by applicable law or agreed to in writing, software\n",
        "distributed under the License is distributed on an \"AS IS\" BASIS,\n",
        "WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.\n",
        "See the License for the specific language governing permissions and\n",
        "limitations under the License."
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "_9XMWz7ir3dX",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "import numpy as np\n",
        "import tensorflow as tf\n",
        "import matplotlib.pyplot as plt\n",
        "\n",
        "# NOTE: Due to the constraints of the /google_research/ repo, this package is\n",
        "# assumes that it is being called from the folder above this one. Because\n",
        "# of this, this notebook needs modify the system path to include the path\n",
        "# above it. This should not be necessary when using this library elsewhere\n",
        "import os\n",
        "import sys\n",
        "module_path = os.path.abspath(os.path.join('..'))\n",
        "if module_path not in sys.path:\n",
        "  sys.path.append(module_path)\n",
        "import robust_loss.general\n",
        "import robust_loss.adaptive\n",
        "\n",
        "# Construct some regression data with some extreme outliers.\n",
        "np.random.seed(1)\n",
        "n = 50\n",
        "scale_true = 0.7\n",
        "shift_true = 0.15\n",
        "x = np.random.uniform(size=n)\n",
        "y = scale_true * x + shift_true\n",
        "y += np.random.normal(scale=0.025, size=n)\n",
        "flip_mask = np.random.uniform(size=n) > 0.9\n",
        "y = np.where(flip_mask, 0.05 + 0.4 * (1. - np.sign(y - 0.5)), y)\n",
        "\n",
        "x = tf.convert_to_tensor(x, tf.float32)\n",
        "y = tf.convert_to_tensor(y, tf.float32)\n",
        "\n",
        "\n",
        "class RegressionModel(tf.Module):\n",
        "  # A simple linear regression module.\n",
        "  def __init__(self):\n",
        "    self.w = tf.Variable(0.)\n",
        "    self.b = tf.Variable(0.)\n",
        "\n",
        "  def __call__(self, z):\n",
        "    return self.w * z + self.b\n",
        "\n",
        "\n",
        "def plot_regression(regression):\n",
        "  # A helper function for plotting a regression module.\n",
        "  x_plot = np.float32(np.linspace(0, 1, 100))\n",
        "  y_plot = regression(tf.convert_to_tensor(x_plot)).numpy()\n",
        "  y_plot_true = x_plot * scale_true + shift_true\n",
        "\n",
        "  plt.figure(0, figsize=(4, 4))\n",
        "  plt.scatter(x, y)\n",
        "  plt.plot(x_plot, y_plot_true, color='k')\n",
        "  plt.plot(x_plot, y_plot, color='r')"
      ],
      "execution_count": 0,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "WrFNwRlbr3da",
        "colab_type": "code",
        "outputId": "33ce112e-cbb4-42c2-ce67-0b371d910c29",
        "colab": {}
      },
      "source": [
        "# Fit a linear regression using mean squared error.\n",
        "regression = RegressionModel()\n",
        "variables = regression.trainable_variables\n",
        "optimizer = tf.keras.optimizers.Adam(\n",
        "    learning_rate=0.01, beta_1=0.5, beta_2=0.9, epsilon=1e-08)\n",
        "\n",
        "for epoch in range(1001):\n",
        "\n",
        "  def lossfun():\n",
        "    # Hijacking the general loss to compute MSE.\n",
        "    return tf.reduce_mean(\n",
        "        robust_loss.general.lossfun(y - regression(x), alpha=2., scale=0.1))\n",
        "\n",
        "  optimizer.minimize(lossfun, variables)\n",
        "  if np.mod(epoch, 50) == 0:\n",
        "    print('{:<4}: loss={:0.5f}'.format(epoch, lossfun()))\n",
        "\n",
        "# It doesn't fit well.\n",
        "plot_regression(regression)"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "0   : loss=13.92873\n",
            "50  : loss=1.78875\n",
            "100 : loss=1.67389\n",
            "150 : loss=1.67370\n",
            "200 : loss=1.67375\n",
            "250 : loss=1.67372\n",
            "300 : loss=1.67424\n",
            "350 : loss=1.67373\n",
            "400 : loss=1.67391\n",
            "450 : loss=1.67377\n",
            "500 : loss=1.67382\n",
            "550 : loss=1.67370\n",
            "600 : loss=1.67370\n",
            "650 : loss=1.67370\n",
            "700 : loss=1.67521\n",
            "750 : loss=1.67386\n",
            "800 : loss=1.67425\n",
            "850 : loss=1.67382\n",
            "900 : loss=1.67464\n",
            "950 : loss=1.67372\n",
            "1000: loss=1.67403\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQgAAAD8CAYAAACLgjpEAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xt8zuX/wPHXtZmZ4xyGzCEl50gWSSWn8HWMIZToQKJCRCVnIb8OzkLIqZXTspzPohiSM4Vy2JINM2zsdP3++Gyz3btv7tm9+7D7/Xw89njYfV/73O9t7veuz3V4X0prjRBCmOPh6ACEEM5LEoQQwiJJEEIIiyRBCCEskgQhhLBIEoQQwiJJEEIIiyRBCCEskgQhhLAol6NeuFixYvrhhx921MsL4bYOHDgQqbX2s6atwxLEww8/zP79+x318kK4LaXUOWvbyi2GEMIiSRBCCIskQQghLJIEIYSwSBKEEMIiSRBCCIskQQghLJIEIYSwSBKEEMIih62kdITgg2FM2nCK8KhYSvn6MLhZJdrV8nfa62YHV4pV3Js9fpdukyCCD4bx0cojxMYnAhAWFctHK48AZOmHml3XzQ6uFKu4N9Pf5dnjf/DhlatAPZv+Lt3mFmPShlOpP8wUsfGJTNpwyimvmx1cKVZxb2l/lzePbuXS0iH8u2muzX+XbtODCI+KzdTjjr5udnClWMW9hUfFonUSUTsXEb1nGd5lH8e3wWs2/126TYIo5etDmJkfXilfH6e8bnZwpViFwdI4Q4m8cHjJeGL//I38NZtRpOnbKE8vm/8u3eYWY3CzSvh4eaZ7zMfLk8HNKjnldbODK8Uq7o4zhEXFork7ZjR33T4uLx1C7F97KdzoLYo064fy9MqW36Xb9CBSBm5sPeqbXdfNDq4UqzA/ZhR17jjvTB5HHuIZNmUBW2/6Z+vvUjnqbM6AgAAtBWOEsKz80DWkfXfeOr6DK+sm45mvMF/MXcp7gY0e6LpKqQNa6wBr2rrNLYYQriZlPEHrJKJ2LSEyZBK5S1agZPcvmXkojuCDYdkegyQIIZzU4GaV8CaeyNWTuL77e/JVb0KJzuPwzFvIbtPTbjMGIYSrqVNCkRQykpiTh/F9oScF67RHKZX6vD2mp63qQSilmiulTimlTiulhpp5vpBSKkQpdUgpdUwp1dP2oQrhPg4cOMBTTz3Fv+dOU/nV0RSq2yFdcgD7TE/ftwehlPIEpgNNgYvAPqXUaq318TTN+gLHtdatlVJ+wCml1BKtdVy2RC2EizBdx9Cwsh/bTkbcc+ZhxYoVvPrqq/j5+fHrr79yNrFoumXVYL/paWtuMeoAp7XWZwGUUkFAWyBtgtBAAWWkuPzAVSDBxrEK4VLM7X1ZvOd86vOme2G01nz22WcMGzaMevXqsWrVKkqUKEGN5PaOmJ62JkH4AxfSfH4RqGvSZhqwGggHCgCdtdZJphdSSvUCegGULVv2QeIVwmWYW8dgKmWwsXmVorz55pssWbKEbt26MXfuXPLkyZParl0tf4esV7FmDEKZecx08UQz4A+gFPAEME0pVTDDF2k9W2sdoLUO8POz6mAfIVyWtYOIF8LCadiwIUuWLGHs2LEsWrQoXXJwJGt6EBeBMmk+L43RU0irJzBBG6uuTiul/gYqA6E2iVIIF+Sb14trMfH3bBN3+W+urBxDZNwNli9fTocOHewUnXWs6UHsAx5TSpVXSuUGXsa4nUjrPNAYQClVAqgEnLVloEK4mvstUo75ay+XFg+mgLcHu3btcrrkAFb0ILTWCUqpfsAGwBOYp7U+ppR6O/n5WcAYYIFS6gjGLckQrXVkNsYthNO7Hmu+96C15kboCq7t+I4KVWuwY+NaSpUqZeforGPVQimt9Vpgrcljs9L8Oxx40bahCeHazG2v1wnxxG6fxbUDG+jUqRMLFizAx8d5t9vLUmshsonp9vrEmOtE/DiMiAMbGDFiBN9//71TJweQpdZCZJu02+v/+eskV1aNRcdc5fvvv+fll192cHTWkQQhRDZqV8sf70uH6fzZUArny8dP63dSp04dR4dlNbnFECKbaK35+uuvadWqFRUqVCA0NNSlkgNIghAiW8TFxdG7d28GDBiAb5VnuPLCJ7y85K/sqeGgNSxcCIcP2/zSkiCEsLGrV6/SvHlz5syZQ9FnOpO/5Yeo3HlS917YNEkcOwYvvACvvQazZ9vuuskkQQhhQydPnqRu3brs3r2bCh2Hkv+5V1Hq7tvMZoVebt6EDz+EJ56Ao0dh7lyYMiXr1zUhCUIIG9m0aRNPP/000dHRbN++nYRHnjXbLkuFXrSGlSuhShWYNMnoOZw6BW+8AR62fztLghDCBqZPn06LFi0oW7YsoaGh1KtXz2JBlwcu9HLmDLRsCR06QJEisHu30XMoViwLkd+bJAghsiAhIYG+ffvSr18/WrRowe7duylXrhxgw3NIbt+G0aOhWjXYtQu++goOHIBnnrHVt2GRrIMQ4gFFRUXRqVMnNm3axKBBg5gwYQKenncTgk3OIdmwAfr1g9On4eWX4YsvwI77NiRBCPEATp8+TatWrTh79izffvstr7/+utl2D1zo5eJFGDgQli2DihVh0yZo0iSLUWeeJAghMmnbtm106NABDw8PNm/ezPPPP2+7i8fHG7MRI0dCQgKMHQuDBoG3t+1eIxNkDEKITJg9ezYvvvgiDz30EKGhobZNDrt2wZNPGgmhQQNjjcMnnzgsOYAkCCGskpiYyIABA+jduzdNmjTh119/5ZFHHrHNxSMioGdPeO45iI6GVasgJARsdf0skAQhxH1ER0fTunVrvv76a9577z1CQkIoVKhQ1i+clATffAOVKsHixTB0KBw/Du3agTJXCtb+ZAxCuCXT8yoszS6cPXuWNm3acOrUKWbNmkXv3r0f6DoZHDgA77wDoaHGUukZM4zFT05GEoRwO+bOq0h7PkWKX375hfbt25OYmMiGDRto1KjRA10nnago+PRTIyH4+Rk9h65dnabHYEpuMYTbMXdehekeiQULFtC4cWMKFy7Mnj17MiQHa6+TSmtYsgQqVzaSwzvvwMmT0K2b0yYHkB6EyCEy09W3tBciPCqWxMREPv74Yz7//HMaN27MsmXLKFy4cKavk87x49C3L2zfDnXqwJo1ULu21d+bI0mCEC4vs119c8VkAUr4QPv27Vm9ejV9+vRh8uTJeHl5WUw+lq6Tutfi1i0YM8ZY/Zg/P8yaBW++CZ6eGb7GWckthnB5merqY36PRK5bkfy3ZDBr1qxh2rRpzJgxIzU5fLTyCGFRsWhIV9PB4l6LFytCcDBUrQoTJ8Irrxg7Lnv3dqnkANKDEC4s5S+7ub/iYPkWwHSPRIHrZzn/42hIjGft2rW8+OLdExzulXx2D22U7jqlfH0YUS0PLw5/27iNqF4dfvkFnjW/7dsVSIIQLsn0tsIc37xeFp9L2SOxZMkSXn99EKVLl+bnn3+mislU4/3GGVL3Wty5Y9Rn6DQOcuWC//s/eO898LIcgyuQWwzhkqw5Ofv2PZ5PSkpi2LBhvPLKK9SrV4/Q0NAMyQEs125I9/jmzVCjhjF92aoVnDgBH3zg8skBJEEIF2VNVabY+CSz9R9v3bpFp06dGDduHG+88QYbN26kaNGiZq9xz5oOYWHGFuymTY1VkevXG7svS5d+sG/KCUmCEC7J2qpMpgOVFy9e5Pnnn2flypV8+eWXzJkzh9y5c1v8+na1/Bnf/nH8fX1QgL+vDxPaVKHd9h+NNQ3BwTBqFBw5As2aZeVbckoyBiFc0uBmle47BgHpexr79u2jbdu23Lhxg5CQEFq2bHnPrzWd3vyq8xO0iz0Hb7YxSsy3aAFTp8Kjj9rke3JGkiCESzKdiQDQZtqlDFT++OOPvPbaa5QoUYLffvuN6tWr3/P6poOgMeGXSOw5EQ5tNG4hVqyAl15y6lWQtqC0NvdjzX4BAQF6//79DnltkfM8MWojUbHxGR4vlCcX7dUeRo4cSf369Vm5ciXFixe/7/XqT9hKWFQsSifR+dBGhuz4jvxxMfz4bAe6rfnWWPjkopRSB7TWAda0lR6EyBGum0kOSfF3OP3TREae/IXu3bsze/ZsvK0svhIeFUvV/84yduN0ngw/xd4y1RnWtA+n/crRzYWTQ2ZJghA5gumy54SbV4lYOYa4S6fp/t7H/FWyEZVHbLZuS3Z0NBN/mUeH34K55lOAAS0HsqpaQ1AK/wctWe+iZBZD5AhppyPj/jvDpYUDiY+8QOCQrwgt+Bzh129nWCqdgdYQFASVK9Pxt1X88OT/aPTWN6yq3giUerCS9S5OEoTIEVKmI/OEHeDSkg/x9PDgi4U/EV7ocev2aZw6Zaxn6NIFSpVChYaSd84sCpT0S53eHN/+8QerUO3C5BZD5Ahaa06sX8ipxSOoW7cuwcHBlCxZkilD15htnzr9GRMD48YZy6Tz5TNqNfTqBZ6etOMehV/chCQI4fJu377NW2+9xeLFi+nSpQvz5s0jT548gOWt3aV8fWD1amO/xLlzxhmXn38OVsxwuBO5xRAu7fLlyzRu3JjFixczduxYlixZkpocwPxS6UdvRbJ83QRo29aYrtyxAxYskORghvQghMs6cuQIrVq1IiIigmXLlhEYGJihTdoFVRFXovngcAhv7FxKLk9Po8fQv3+O2FSVXSRBCJcUEhJC165dKViwIDt37iQgwPK6n3a1/Gl39ST0HWoMRnboYByAW6aMHSN2TVbdYiilmiulTimlTiulhlpo84JS6g+l1DGl1A7bhimEQWvNF198Qdu2balUqRKhoaH3TA78+69RGLZJE+NYu7VrYflySQ7W0lrf8wPwBM4AjwC5gUNAVZM2vsBxoGzy58Xvd93atWtrITLjzp07umfPnhrQgYGB+tatW5Ybx8drPWWK1gULau3trfXw4VrHxNgvWCcG7Nf3eX+mfFjTg6gDnNZan9VaxwFBQFuTNl2BlVrr88lJ53IW85YQ6URGRtK0aVPmz5/P8OHD+eGHH8ibN6/5xnv2wFNPGTMU9eoZW7FHjQIf91oFaQvWJAh/4EKazy8mP5ZWRaCwUmq7UuqAUqq7uQsppXoppfYrpfZHREQ8WMTC7Rw/fpw6deqwd+9eli5dyqhRo/DwMPNf98oVYw1DvXrGeZfLlsG6dfDYY/YPOoewJkGY289qugU0F1AbaAk0Az5VSlXM8EVaz9ZaB2itA/z8/DIdrHA/69evp169esTExLBjxw66dOmSsVFSEsyfbxRwmTfPKPd24gQEBub47djZzZoEcRFIO6JTGgg302a91vqW1joS2AnUtE2Iwh1prZkyZQotW7akfPny7Nu3j7p162ZseOiQcSr2668bh+AePGgUjC1QwP5B50DWJIh9wGNKqfJKqdzAy8BqkzY/Ac8ppXIppfICdYETtg1VuIv4+Hj69OnD+++/T+vWrdm1axdlTGcdoqNh4EDjhKo//zR6Djt3wuOPOyboHOq+6yC01glKqX7ABowZjXla62NKqbeTn5+ltT6hlFoPHAaSgLla66PZGbjIma5evUrHjh3ZunUrH330EWPHjk0/3qC1MbYwYIAxhdmrF3z2GRQp4rigczJrpzts/SHTnMLUyZMn9WOPPaZz586tv/vuu4wNTp3SumlTrUHrJ5/Ueu9e+weZA5CJaU5ZSSmcwubNm+nYsSNeXl5s3bqV+vXr330yNhbGjzeOscuTxygU26ePyx1j54pks5ZwuJkzZ9K8eXNKly5NaGho+uSwZg1Uq2Ycgtuxo7FUul8/SQ52Ij0Ike0snY6dkJDAwIEDmTp1Ki1btmTp0qUULFjQ+KLz5+H9941zJ6pUga1boWFDx34jbkgShMhWpuXjU0q+3bpxnYXjBrBx40Y++OADJk6ciKenJ8TFGRupRo82LjB+vDFbcY/DbUT2kQQhspW5MzSjL1/grcDexF8LZ+7cubzxxhvGE9u3wzvvGIuc2raFyZOhXDn7By1SSYIQ2cr0DM3b5w8TsWo8AFs3b6ZBgwbw338waBAsXgzly0NIiHEIrnA4GaQU2SrtGZo3Dm3kvx8+xSNvIWr2nU6DZ5+F6dONFZA//micjn3smCQHJyIJQmSrwc0qkccTrm6dy9X1U8hTtgblX/+aL2sWhjp1jBmJp54yzrocPVp2XDoZucUQ2arRowXIv/MrbuzbQoHaranX7DVm/fkT5b9YDCVLGudQdOokm6qclCQIkW3+/vtvWrduzcmTJ5kxfTp98uaFD/vA1avGFOaoUZAyrSmckiQIYVMpax7OHj1AZPA4vD3gl5kzqbdwIezaZdRqmDkTaspmX1cgCULYTMqah4iDG7myfiq+BYrxebka1O3TB3x94dtvoUcPMFfsRTglSRDCZj5fd4Lwzd8SvWc5XYqV5YvYmzx0aCM/PdWStuu+g6JFHR2iyCRJEMImbt68yaEFn/LQX3tYXqAYTSPPc6z4I/R96SMO+lehrSQHlyQJQmTZ+fPnCWzVisF/HWGohyfxd24xqvFbLHyyFYkenhTOKwfTuCpJECJL9uzZw5QWLQi6fp1HgJ8r12dUwzeJyH+3gMvN2wkEHwxz+4NwXZGMFokH9tO0aYTXr8/SqCj8H34YtmxhWKdP0iUHgPgkzaQNpxwTpMgSSRAi05Lu3GFDkyY0fvdd/gfc+uQTvE+cgEaNuB4bb/ZrTPdkCNcgtxgiU25v2sTlwECaRUdzsGxZqm3aRJ6Kd084KOXrQ5iZZJB2T4ZwHdKDEICxhqH+hK2UH7qG+hO2EnwwLH2Dy5e51akTeV58kaToaIJ79uSJf/4hd8X0x58MblYJH6/01Z58vDwZ3KxSdn8LIhtID0JYLOoC0K5GSZg9m4QhQ/C6cYMvvLyosmQJ7Tp2NHutlIFIcxWkhOuRBCHMFnWJjU8k5NvVtNs7D/bv5xcPD8Y89BBfr19PjRo17nm9drX8JSHkEJIgRIYBxIK3bzJo5yJeObiWG/nz0Qs4X7cuq4KDKV68uGOCFA4hYxDi7gCi1rQ/uoUtc96m2x/rmFvMn9I3b+L16qts3bZNkoMbkgQhGNysEo9HXeSH7z/iyzVfcaFAUZ4tWobekRcZ+tlnfPfdd3h7ezs6TOEAcovh7m7epF3QFNp8+yU3vHwY9MzLTDu6BeJusnLlSl566SVHRygcSBKEu9LaOHPi/ffhwgU8Xn+dX194gVl9+uBXuDCrV2+gVq1ajo5SOJjcYrijM2egZUto3x58fdG//MLEihVp+dprVKtWjdDQUEkOApAehHu5fZsTA4bxyNwpxHvkYl7LPpT6qD8hs8excOFCOnfuzPz58/GRwrEimSQId7FxIzff7E2VC/8QUvk5xjZ6g3BPL668/BKxF48zcuRIhg8fjpLisSINSRA5XVgYDBgAy5Zxtag/b3caw67ytYiL+IfLK8aQdOsa5QI/ZsSIEY6OVDghSRA5VXw8TJ0KI0ZAQgKMHk3T649zJ5cXMWf2Ebn6czxy+1Ci6wR4qKLUaxBmySBlTrRrF9SuDR98AA0aGKdVffopRYsWIDp0FRHLR+NVuBQlX/0S74eMzVZSr0GYIwkiJ4mIgNdfh+eeg+vXjWnMkBB45BHi4uLIt+9brm37lrwV61Gi60RyFSyW+qVSr0GYI7cYOUFSEsydC0OHwo0bMGSIcc5lvnwEHwzjs5WhHP5uOHcuHKVI/c7kr98NpdL/bfCVupHCDEkQru733+Gdd2DvXuN2YsYMqFoVMLZxD5y9hgtBI0m4EUmx1oPIV/UFs5fR2o4xC5chtxiuKioK3n3XOPj2n39g0SLYti01OQB8PHUx/8wbQFJ8LCW7jLeYHACLpeKEe5MehJNLOcoutfjKixVpd3Int9/rj9fVKyyp1YIlrXrRp1pt2iWvYdBaM23aNE4s+Bgvv3IU7/ApuQreeyemlIQT5liVIJRSzYHJgCcwV2s9wUK7p4A9QGet9XKbRemmTCs95Tn9JyVnDYRzh/izVEU+7v4JR0tWgDukVoBqWb0477//PjNnzqRwlWfI33wAHrnTv/kVkPaOQkrCCUvumyCUUp7AdKApcBHYp5RarbU+bqbdRGBDdgTqjlIqPfnE3ebd34J4MzSYWC9vhjXry9IaL5Lkcbf2Y2x8IuOD9zNt8FS2bNnChx9+SN1O/fgk+Fi6alE+Xp50qO3PtpMRUhJO3Jc1PYg6wGmt9VkApVQQ0BY4btLuXWAF8JRNI3Rj4ddiaHp6LyM2f0Pp6AiWV2/M+Bd6ciWfb4a28VfD+H3FaNSNy8yfP58ePXoA4OHhIfUhxQOzJkH4AxfSfH4RqJu2gVLKH3gJaIQkCNv4+28WrR7Hsyf3cLJYOTp2ncC+MtUB8FSKxDTTDrHnDhEZPB4PT0+2bNnCc889l/qc1IcUWWHNLIa53Tumk2JfA0O01olm2t69kFK9lFL7lVL7IyIirI3Rvdy5A+PGQdWqPH3uMBObvEmrHpNTk4OPlydd6pZJLS1/4491XP5xOLkKFGVa0Lp0yUGIrLKmB3ERKJPm89JAuEmbACAoeSdgMeB/SqkErXVw2kZa69nAbICAgACZeTe1eTP07Qt//gmBgeT66isqRShKmLlFqOVfkHf7D+TqryvxrVSX6XMW0PW5yo7+DkROo7W+5wdGEjkLlAdyA4eAavdovwAIvN91a9eurUWysDCtO3fWGrR+9FGt1627Z/OoqCjdvHlzDej+/fvrhIQEOwUqcgJgv77P+zPl4749CK11glKqH8bshCcwT2t9TCn1dvLzs2yfttxEQgJMn24si46Lg5EjWd3sFSZuP0f49jVmBxXPnDlD69at+euvv5g9ezZvvfWWA78BkdNZtQ5Ca70WWGvymNnEoLXukfWw3MBvv0GfPnDoEDRrBtOmEXzDx/IJV7X82blzJ+3btycpKYmNGzfSsGFDR34Hwg3IUmt7i4yEN9+EZ54x/r18OaxbBxUqWDzhatKGU8ybN48mTZrg5+dHaGioJAdhF7LU2l6SkmDePGOnZXQ0DB4Mw4dD/vypTcxtudZJiRxdNZ03QldSqEJthn2zkAoVKtgzcuHGJEHYwx9/GLcTe/YYtRpmzIDq1TM0K+XrQ1iaJJF0J4bIkEnEntlHgSdbUqhxL8ZuOk++AoVkbYOwC7nFyE7R0dC/v1Hd6cwZWLAAduzIkByCD4ZRf8LWdMkh4fplLi35kNizByjStA9FmvZBeXim3nIIYQ/Sg8gOWsMPP8DAgXDpEvTuDZ99BoULZ2hquiEL4PbFE0SsGodOjKd4x1H4lE9/RoVUfxL2IgnC1k6dMhY7bdli9Bx++smo2WCB6cDkzWPbuLJuMrkK+lG8w3i8ipbJ8DWyNVvYiyQIW4mJMXoJn38OefMa6xt69wZPz3t+WUpvQOskonYuInrPMrzL1sCv3Ud4+hTI0F62Zgt7kgRhCyEh8N57RmWnV1+FSZOgRAmrvrSUrw8XLl8jcs0XxP75G/lrNjPGGzyNX41S4OvjRVRMvOzGFHYnCSIr/vnHOPx29Wqj1Nv27UZdyEzo+UQB+vXsxZ3//qZwo7coENAm/elWGg4Of9GmYQthLUkQDyIuDr74AsaMMf7ET5xonF7llbnK0KGhoXzasy0eN25SMnAE3o/UztBGxhuEI8k0Z2Zt2wY1a8LHH0Pz5nDiBHz4YaaTQ1BQEA0aNMDHx4d9e/cw86M3Urdwp5DxBuFokiCsdekSvPIKNGoEcXH8NuU76tfpR/kZR6g/YSvBB8OsukxSUhIjRoygS5cuBAQEsHfvXqpVq0a7Wv6Mb/84/r4+KMDf14fx7R+X8QbhUHKLcT8JCTBzJgwbBrdvw/DhrG7+KkPWniY23piBMN1UZUlMTAw9evRg2bJl9OjRg1mzZuHt7Z36vFR/Es5GehD3sncv1KljzFA8/TQcPQqjRjFxx3mLm6osCQ8Pp0GDBixfvpyJEycyb968dMlBCGckPQhzrl6Fjz6COXPgoYfgxx8hMNAYkMTySkZLj//++++0adOGqKgogoODadOmTbaFLoQtSQ8iraQkmD8fKlWCb7819lGcPAkdO6YmB7A8s2Du8RUrVvDss8/i4eHB7t27JTkIlyIJIsXhw/D888bp2JUqGWdefvklFMi4mnFws0r3nXHQWjN27FgCAwOpWbMmoaGh1KxZM9u/DSFsSW4xbtyAkSNh8mRjM9X8+dC9O3hYzp0pA4lpz5toWNmPSRtOMeCHPyiZz5Pcv33DznXBdOvWjblz55InTx47fUNC2I77JgitjWpO/fvDv/9Cr17GXooiRaz68rQzDml3ZCbeusbvi8YSF36Kru98yKJpE9KvjBTChbhngvjrL+jXDzZuhFq1YMUKY5biAaXsyIy7fJbLy8eQFBuNX7uP+adMY0kOwqW5V4KIjYXx442l0XnywJQpRqWnXFn7MYRHxRLz114iQybh4Z2PEt0m4l2ygtRtEC7PfRLE2rVGr+Hvv6FrV2MvRcmSWb6s1hp96Cci1s8ld8kK+LUfRq4CRQHZRyFcX85PEOfPG+MMq1ZB5cqwdSvYqCL0nTt36NOnD+fWz6dA1efxbf4eHl7GYKTsoxA5Qc6d5oyLM24lqlSB9euNW4tDh2yWHCIiImjatCnz589nxIgRfLdoCWX8Css+CpGj5MwexI4d8M47cPw4tG1rTGGWK2ezyx89epTWrVtz6dIlgoKC6Ny5MwAvPVnaZq8hhDPIWT2I//4z1jC88IJRAi4kBIKDbZoc1q5dyzPPPMPt27fZsWNHanIQIifKGQkiMdE4a6JSJQgKgk8+gWPHoFWrB7pcShn68kPXpG7l1lrz1Vdf0bp1ax599FFCQ0OpU6eOjb8RIZyLMg77tb+AgAC9f//+rF9o3z5jqvLAAWjc2CgWW+nBBwfNlaHP46EpdyqIjSuX0L59exYuXEi+fPmyHrsQDqCUOqC1DrCmreuOQVy7ZlR1+uYbY7oyKAg6dUq3qepBmJahT4y9wT/Bn3Hq/BE69OxHeIW2VB+zXQrICrfgeglCa1i40Djb8soVo2jsqFFQsKBNLp92cVP8lQtcXj6ahBsRFG31Acf8mxAbfQewvkiMEK6pOLjzAAAKHUlEQVTMtcYgjh41qkb36AEVKhi3FV99ZbPkAHcXN8X+fZB/Fw0iKS6Wkl3GU6h6o0wXiRHC1blGgtDaOBX7iSeMwcc5c2DXLuNzGwo+GEZMXAI3fv+Zy8tGkKugHw91/xLfh6uTaGGsRpZTi5zMNW4xlIKbN6FnT2PBU7FiWbpc8MGwdFu1U1Y8Dl3+B+HrZ3Lj9zX4VKhDsVaDKOJbiJFtqjFpw6l0h+umkOXUIidzjQQBMG1algcgIeMsRcpYQq6EGM4FjeX2PwcpWKc9vg1eQ3l4ks87V+oYg+nshiynFjmda9xigE2SA2ScpQCI/u88J2e/z+3zRyja4n0KN3wd5WFUjAqLiqX+hK0AUpZeuB3X6UHYiOmYwe1zh4kI/gyUByVeHkueMtUzfE1KL2N8+8fZPbSRvUIVwuFcpweRRSmrI9MONd74Yz3//fgpnvkK8/DrX1H4Ecs1I2XGQrgjt0gQKeMOKYOMOimRq1vmcHXDNPKUe4KSr/4fHoVK0qG2P/73GHSUGQvhbtwiQaQdd0i6E8PlFaO5sf8nCtRuQ/HA4Xh45yM+UbPtZAS7hzaymCRkxkK4G6sShFKquVLqlFLqtFJqqJnnuymlDid//KqUcqr67il/+eOjLnFp0SBu//MHRZr1pUiTXqmDkWnbWVPWXgh3cN9BSqWUJzAdaApcBPYppVZrrY+nafY30EBrfU0p1QKYDdTNjoDNMbeuIe3sQilfH84c2UfEqs9AJ1G802h8ymXMYSk9BHNl7WXfhXBH1sxi1AFOa63PAiilgoC2QGqC0Fr/mqb9HsBulVMsrWuAu2/02ncO8WvQMHL5lqB4h+F4Fcn4RjftIchBukJYd4vhD1xI8/nF5McseQNYl5WgMsPcuoaUGYfExESGDBnC1JEDqRHwNE/2m07uIsZA5CtPl5U1DULchzU9CHMrlMxuTFBKNcRIEM9aeL4X0AugbNmyVoZ4b5ZmFi5evkr79u1ZvXo1ffr0YfLkyXh5ednkNYVwF9b0IC4CZdJ8XhoIN22klKoBzAXaaq2vmLuQ1nq21jpAax3g5+f3IPFmYG5mISH6MpHfD+Hnn39m6tSpzJgxQ5KDEA/AmgSxD3hMKVVeKZUbeBlYnbaBUqossBJ4VWv9p+3DtMx0xuFO2AkuLRwINyNZu3Yt/fr1s2c4QuQo973F0FonKKX6ARsAT2Ce1vqYUurt5OdnAcOBosCM5KPmEqwtafWg0s5cFPLxIo+XBxf3beTKuimUeKgUWzeuo0qVKtkZghA5nkvWpDSdudA6iVu7l3JldxANGjRgxYoVFC1a1JbhCpFj5PialOlWRsbd5sqaL4n581eKP/U/Nm5cRe7cuR0coRA5g0smiJSZi4ToSCJWjiHuv7MUbvQmPgFtJTkIYUMumSBK+fpw9sQhIlaOJSkuFr/A4eR99Kl7brQSQmSeS27Wekb9yX9Lh4KnFyVfmUTeR5+SvRJCZAOX6kForRk9ejT/N3IkVZ6oQ/6WQ4hI8Ja9EkJkE5dJELGxsfTs2ZMffviB7t27M3v2bLy9vR0dlhA5mkskCK01LVq0YOfOnUycOJHBgwejbFSjUghhmUskCKUUffv2pX///rRr187R4QjhNlxikDL4YBhfnynKgD1eqadtCyGyn9OvpDR32raXhyJ/nlxExcTLAKUQmZSjVlKaq/cQn6S5FhMPyCG6QmQnp7/FsKaStJSkFyJ7OH2CsLaStJSkF8L2nD5BmKswbY6UpBfC9px+DMK0wnQhHy9uxSUQn3h3cFWWWQuRPZw+QUDGCtP3K3MvhLANl0gQpqQkvRD24fRjEEIIx5EEIYSwSBKEEMIiSRBCCIskQQghLJIEIYSwSBKEEMIiSRBCCIskQQghLJIEIYSwSBKEEMIil9yLIYQtOHrTn6Nf3xqSIIRbMq11au/ShY5+fWvJLYZwS+ZqndqzdKGjX99akiCEW7JUotBepQsd/frWkgQh3JKlEoX2Kl3o6Ne3liQI4ZbM1Tq1Z+lCR7++tWSQUrgl01qn9p5FcPTrW8vpT9YSQthWZk7WklsMIYRFkiCEEBZJghBCWGRVglBKNVdKnVJKnVZKDTXzvFJKTUl+/rBS6knbhyqEsLf7JgillCcwHWgBVAW6KKWqmjRrATyW/NELmGnjOIUQDmBND6IOcFprfVZrHQcEAW1N2rQFFmrDHsBXKfWQjWMVQtiZNQnCH7iQ5vOLyY9lto0QwsVYkyCUmcdMF09Y0walVC+l1H6l1P6IiAhr4hNCOJA1KykvAmXSfF4aCH+ANmitZwOzAZRSEUqpc5mKFooBkZn8GnuS+LJG4ssaa+MrZ+0FrUkQ+4DHlFLlgTDgZaCrSZvVQD+lVBBQF7iutf73XhfVWvtZG2QKpdR+a1eAOYLElzUSX9ZkR3z3TRBa6wSlVD9gA+AJzNNaH1NKvZ38/CxgLfA/4DQQA/S0ZZBCCMewarOW1notRhJI+9isNP/WQF/bhiaEcDRXW0k529EB3IfElzUSX9bYPD6H7eYUQjg/V+tBCCHsyCkThLPv/bAivm7JcR1WSv2qlKrpTPGlafeUUipRKRXobPEppV5QSv2hlDqmlNrhLLEppQoppUKUUoeSY7PrgLxSap5S6rJS6qiF52373tBaO9UHxkzJGeARIDdwCKhq0uZ/wDqMBVpPA3udLL5ngMLJ/27hbPGlabcVY/A50JniA3yB40DZ5M+LO1FsHwMTk//tB1wFctvx5/c88CRw1MLzNn1vOGMPwtn3ftw3Pq31r1rra8mf7sFYOGYv1vz8AN4FVgCX7RgbWBdfV2Cl1vo8gNbaXjFaE5sGCiilFJAfI0Ek2Ck+tNY7k1/TEpu+N5wxQTj73o/MvvYbGBndXu4bn1LKH3gJmIX9WfPzqwgUVkptV0odUEp1d6LYpgFVMFYKHwHe11on2Sc8q9j0veGMRWtttvcjm1j92kqphhgJ4tlsjcjkZc08Zhrf18AQrXWi8YfQrqyJLxdQG2gM+AC/KaX2aK3/dILYmgF/AI2AR4FNSqlftNbR2RybtWz63nDGBGGzvR/ZxKrXVkrVAOYCLbTWV+wUG1gXXwAQlJwcigH/U0olaK2DnSS+i0Ck1voWcEsptROoCWR3grAmtp7ABG3c8J9WSv0NVAZCszk2a9n2vWGvwZVMDMLkAs4C5bk7UFTNpE1L0g/EhDpZfGUxlp0/44w/P5P2C7DvIKU1P78qwJbktnmBo0B1J4ltJjAy+d8lMPYnFbPz7/hhLA9S2vS94XQ9CO3kez+sjG84UBSYkfxXOkHbaZOPlfE5jDXxaa1PKKXWA4eBJGCu1trstJ69YwPGAAuUUkcw3oRDtNZ22+GplPoeeAEoppS6CIwAvNLEZ9P3hqykFEJY5IyzGEIIJyEJQghhkSQIIYRFkiCEEBZJghBCWCQJQghhkSQIIYRFkiCEEBb9P9j1geIQqkiFAAAAAElFTkSuQmCC\n",
            "text/plain": [
              "<Figure size 288x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "JxKjAb2Cr3dd",
        "colab_type": "code",
        "outputId": "20fad39a-9d49-495f-f7cc-68284b5e5e64",
        "colab": {}
      },
      "source": [
        "# Fit a linear regression, and the parameters of an adaptive loss.\n",
        "regression = RegressionModel()\n",
        "adaptive_lossfun = (\n",
        "    robust_loss.adaptive.AdaptiveLossFunction(\n",
        "        num_channels=1, float_dtype=np.float32))\n",
        "variables = (\n",
        "    list(regression.trainable_variables) +\n",
        "    list(adaptive_lossfun.trainable_variables))\n",
        "optimizer = tf.keras.optimizers.Adam(\n",
        "    learning_rate=0.01, beta_1=0.5, beta_2=0.9, epsilon=1e-08)\n",
        "\n",
        "for epoch in range(1001):\n",
        "\n",
        "  def lossfun():\n",
        "    # Stealthily unsqueeze to an (n,1) matrix, and then compute the loss.\n",
        "    # A matrix with this shape corresponds to a loss where there's one shape\n",
        "    # and scale parameter per dimension (and there's only one dimension for\n",
        "    # this data).\n",
        "    return tf.reduce_mean(adaptive_lossfun((y - regression(x))[:, None]))\n",
        "\n",
        "  optimizer.minimize(lossfun, variables)\n",
        "  if np.mod(epoch, 50) == 0:\n",
        "    loss = lossfun()\n",
        "    alpha = adaptive_lossfun.alpha()[0, 0]\n",
        "    scale = adaptive_lossfun.scale()[0, 0]\n",
        "    print('{:<4}: loss={:+0.5f}  alpha={:0.5f}  scale={:0.5f}'.format(\n",
        "        epoch, loss, alpha, scale))\n",
        "\n",
        "# It fits!\n",
        "plot_regression(regression)"
      ],
      "execution_count": 0,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "0   : loss=+1.30600  alpha=1.00500  scale=0.99369\n",
            "50  : loss=+0.81370  alpha=1.24806  scale=0.69965\n",
            "100 : loss=+0.41574  alpha=1.46083  scale=0.47681\n",
            "150 : loss=+0.03641  alpha=1.62745  scale=0.31602\n",
            "200 : loss=-0.25180  alpha=1.73791  scale=0.20793\n",
            "250 : loss=-0.37724  alpha=1.64325  scale=0.14308\n",
            "300 : loss=-0.45471  alpha=1.44397  scale=0.11082\n",
            "350 : loss=-0.55913  alpha=1.20640  scale=0.08612\n",
            "400 : loss=-0.71528  alpha=0.94550  scale=0.06072\n",
            "450 : loss=-0.91360  alpha=0.69710  scale=0.04031\n",
            "500 : loss=-1.11517  alpha=0.49126  scale=0.02743\n",
            "550 : loss=-1.27192  alpha=0.33634  scale=0.02091\n",
            "600 : loss=-1.37283  alpha=0.22587  scale=0.01729\n",
            "650 : loss=-1.43477  alpha=0.14965  scale=0.01529\n",
            "700 : loss=-1.47247  alpha=0.09828  scale=0.01377\n",
            "750 : loss=-1.49214  alpha=0.06422  scale=0.01310\n",
            "800 : loss=-1.50190  alpha=0.04189  scale=0.01264\n",
            "850 : loss=-1.50135  alpha=0.02733  scale=0.01239\n",
            "900 : loss=-1.50108  alpha=0.01790  scale=0.01227\n",
            "950 : loss=-1.51328  alpha=0.01181  scale=0.01207\n",
            "1000: loss=-1.51522  alpha=0.00789  scale=0.01212\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAAQgAAAD8CAYAAACLgjpEAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDMuMC4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvnQurowAAIABJREFUeJzt3Xl4jPf6x/H3nZFIUBKkRcRyaq+l2lQQS6m1drVWUdWiRc+vi9Ll0KNaHKeqPYpDqtQh1FrriSW1VTUoLaXUUks4qIYUQZbv74+EJmMmRkwyM5n7dV29rszMd565K9d88jzf7RFjDEopZYuPqwtQSrkvDQillF0aEEopuzQglFJ2aUAopezSgFBK2aUBoZSySwNCKWWXBoRSyq58rvrg4sWLm3Llyrnq45XyWrt27frNGBPsSFuXBUS5cuXYuXOnqz5eKa8lIscdbauXGEopuzQglFJ2aUAopezSgFBK2aUBoZSySwNCKWWXBoRSyi4NCKXygthYuHjR6YfVgFDK00VFQaNG8MYbTj+0y2ZSusKy3XFMiD7I6YuJlAoMYFjLynSsHeK2x80JnlSrytqy70/x2/B3eH79bHaXr8npni/Txsmf4TUBsWx3HG8u2UtiUgoAcRcTeXPJXoB7+oLk1HFzgifVqrK2/Luj+Dz/As/vi2Fx9aa8UrMl/qt/JikwyKm/S6+5xJgQffDWF+OmxKQUJkQfdMvj5gRPqlVl4cIFyvToSPt9MUxo2JuBZWpxcsHbnFkX6fTfpdcExOmLiXf1vKuPmxM8qVZlx6FDULcuVU/+zJB2rzMm6RoXVn9E/pCqBDbu6/TfpddcYpQKDCDOxj9eqcAAtzxuTvCkWlWajH1GbX4/yMQFo/HL78fgfh8wb/MSEg99S6FaLSnafBBi8XX679JrziCGtaxMgK8l03MBvhaGtazslsfNCZ5Uq/qzzyjuYiKd925g4mfDOeFXmAXjphOzZRaJv3xHUNMXKNpyCGLxzZHfpdecQdzsuHF2D35OHTcneFKtKu33dO1GEq9v+Q9Dvv2SrWVr0b9uV479X3/8SeKdT2YRczkkR3+X4qp7c4aFhRndMEYp+6q8toQJqyfR7uctRNVswauh1TkbPRlLwSA+jJzHy12aZuu4IrLLGBPmSFuvucRQyqOcO8eihe/Q5uetfND4WV4sVJT/rZqIX4kKlOgzkak/3GDZ7rgcL0MDQil3s38/hIdT5fwxXu44jPfPHuHStvkUrN6MB7q/j6VAkVwbntaAUMqdrF8P9evDtWvEL15M9KForv68lcDH+1Hsyb8i+XxvNc2N4WmHAkJEWonIQRE5LCIjbLxeRERWiMgPIvKTiPRzfqlK5XEzZkCrVhAayo+RkTz8wgucOX6YKr1HUyT8KUQkU/PcGJ6+4yiGiFiAT4HmwClgh4gsN8bsz9BsMLDfGNNORIKBgyIy1xhzI0eqVspDWK99aVIlmK9/Pp955KFWSRgxAiZMgJYt+apXL3p27UpwcDDbtm3jaEqxTFPkIfeGpx0Z5qwDHDbGHAUQkflAByBjQBjgPkmLuELA70Cyk2tVyqPYWvvyn+0nbr0edzGR0fN3UOe1aZT6OhozaBBjS5bk7T59qFevHkuXLuWBBx6gZnp7VwxPOxIQIcDJDI9PAeFWbSYDy4HTwH1Ad2NMqvWBRGQAMACgTJky2alXKY9ha+1LRsGXfydy8XuUOHuYpH/8g3579jB32jR69epFZGQk/v7+t9p2rB3ikvkqjvRBiI3nrCdPtAT2AKWAh4HJIlL4tjcZM90YE2aMCQsOdujGPkp5rKw6EaucO8ayL16j4oUTvPDk/9FoyRLmzpvHmDFjmDNnTqZwcCVHziBOAaEZHpcm7Uwho37AOJM26+qwiBwDqgCxTqlSKQ8UWMCX+KtJtz3/+JGd/Gv5eK74BdCp1V/ZuHkWlht/sGjRIp566ikXVGqfI2cQO4CKIlJeRPyAHqRdTmR0AngCQEQeACoDR51ZqFKextYk5We+X8Vni0dzIrAkLRv0Zt1/P+G+/D5s3brV7cIBHDiDMMYki8gQIBqwADONMT+JyKD016cB7wGzRGQvaZckw40xv+Vg3Uq5vUuJf549+KSm8E7MZzy3aznrHqzDcyUrcir6EypUq8mmtaspVaqUCyu1z6HFWsaY1cBqq+emZfj5NNDCuaUp5dluLq8vcCORT5b/g2ZHdvDZI215mxuc3TqXbt26MWvWLAIC3He5vc6kVCqHDGtZmXKJv7Nw7nCaHN3FO437MvT8Mc5+v5ZRo0YRFRXl1uEAXrTcW6nc1pFztIwaRmrCH/RpPojFO5Zirv5OVFQUPXr0cHV5DtEzCKVywvLl0KABAf5+fD9xAiu+nUNQfsOWzZs9JhxAA0Ip5zIGPvoIOnbEPPQQ059/nseHDqVChQrExsZSp04dV1d4VzQglHKW5GQYPBhefZXUjh0ZWr06A0eNIrBqfS48/jY95v6SK3s4OJP2QSjlDAkJ0K0bREeTOHQobffuJWbpUorV707BBr0Q8fHI+5DoGYRS9+r4cYiIgA0bODN6NDXXrGHrtm1U6DqCQg17I/Ln18zT7kOiAaHUvYiNhfBwOHmSXe+/T9UPPyQhIYGNGzeS/JcGNt/iSfch0YBQKrsWL4bGjaFAAeYNHUr4W29RpkwZYmNjqVevnt0NXTzpPiQaEErdLWNg/Hjo0gXz8MMMb9yYXmPG0Lp1a7755hvKli0L5I37kGgnpVJ348YNeOkl+OwzbnTuTKeLF1k9axavv/4648aNw2L5MxDywn1INCCUclR8PHTpAjEx/P7SS0SsX8+RY8f47LPPeO6552y+xVUbvTiLBoRSjjhyBNq0gaNH2T9iBA3+/W98fHxYv349jRo1cnV1OUb7IJS6k23boG5dOH+e5S+/TK1//pOSJUsSGxubp8MBNCCUylpUFDRtigkM5P127ejw4Yc0a9aMbdu28Ze//MXV1eU4DQilbDEG3nsPnn6a5LAwepQtyzuzZ/Pyyy+zYsUKihQp4uoKc4X2QSivZH2/ikyjC9evwwsvwJw5/NGxI40OHmTfL78wbdo0Bg4c6Phx8gANCOV1bN2v4tYaiTL+0KkTbNnCr88/z2NLl5KSmkp0dDRNmzZ1/Dh5JCT0EkN5HVv3q0hMSmH+f9andUbGxrJp0CAqzZ5NUNGibN++/bZwyOo4nrTW4k70DELlCXdzqm9rLUT4ib1MXfoBplB+pnbpwuBp03jiiSdYuHAhQUFBDh8nq+c9kZ5BKI9381Q/7mIihj9P9e3tvWC9FqLzvg3MWfA3LhUqwsBatRg8dy4vvvgia9asISgoiGW744gYF0P5EauIGBdz67h5Ya3FnWhAKI93t6f6N9dIiEnl1c1zmLjqI3aUqkT7wgWYuXEjkydPZsqUKfj6+mYZPnlhrcWd6CWG8lg3Lyvi7vJUv2PtEHyuX6PAwOdp9uNGllSpy4Czv5D8RzKrV6+mRYs/7+CQVfh8M6LprTY6iqGUG7EeQbAlsICv7RfOnaP9q73hx2/Z3b07PZcsoXRoKCtXrqRq1aqZmt6pn8HT11rciV5iKI90pztnA1yz9fr+/VC3LmbPHuY99RSPLFhAvfr1iY2NvS0cwDv6GbKiAaE8kiMjBYlJqZk7Ktevh/r1Sb16lRF169Jr8WL69+/P2rVrKVasmM1jeEM/Q1Y0IJRHcvQv+K2OyhkzoFUrkkqUoF1wMBM2bmTixInMmDEDPz8/u+/vWDuEsZ1rEBIYgAAhgQGM7VwjT19WZKR9EMojDWtZ+Y59EABn4q/AG2/AhAlcrFePOkePcubKFVasWEGbNm2yfK/13IqPuj/sNcFwkwaE8kjWuzUBGKs2/knX+HTNR3DgGw43b87DmzdTvEQJvl2/nurVq2d5fG+YRu0IvcRQHqtj7RC+GdGUY+PaUCQg84hF8OXfWTDvTZoc2MZ/W7ak4rp1PBwWRmxs7B3DAbxjGrUjNCBUnnApMenWz1XOHWPZF69R4cJJuoVUpXV0NH369GHDhg3cf//9Dh3PG6ZRO0IDQuUJNzstHz+yk0Vz38CSmkyTIsEsPn2APi+/xS9V+lBl1PpMU6UdOZ6jz+dVGhAqTxjWsjLP/bCGzxaP5lihYtQBdl46T5fhHxFbuCGnL11zaJ1GxuN58/DmTRoQyvOlpNBx9gRG/vdTNpWqSP2Ec5zN58uHX3zF6SI1stWX4O3DmzfpKIbybJcvQ8+esHIlsRERNP/mGx4LD2fZsmWUKFGCT0assvk2R/oS8vo0akfoGYTyXKdOQcOGmDVrmPXYY4R/8w3de/Zk48aNlChRAtC+hHulAaE80/ffQ3g4qYcP81qlSvTbsYMxY8Ywd+5c/P39bzXTvoR7o5cYyvMsXw49e3KjSBHaFirE1l9/ZeHChXTp0uW2pnnh9neupAGhPIcxMGkSvPYa8RUqEBYXx7XAQDZv3kxYWJjdt2lfQvY5dIkhIq1E5KCIHBaREXbaPC4ie0TkJxHZ5NwylddLTobBg+HVVzlUvTqlf/mFoKpViY2NzTIc1L254xmEiFiAT4HmwClgh4gsN8bsz9AmEJgCtDLGnBARx6arKeWIS5egWzdYu5bV1avTdu9enurShdmzZ1OgQAFXV5enOXIGUQc4bIw5aoy5AcwHOli1eRpYYow5AWCMOefcMpXX+vVXiIjAxMTwj4oVabNvH38bOZIFCxZoOOQCR/ogQoCTGR6fAsKt2lQCfEVkI3Af8LEx5gvrA4nIAGAAQJkyZbJTr/ImsbHQvj0pV6/St3hxFp04wbx58+jZs6erK/MajgSE2HjOemVtPuBR4AkgAPhWRLYbYw5lepMx04HpAGFhYdbHUOpPixZB795cDQykcWoqJ41h06ZNhIdb/21SOcmRS4xTQGiGx6WB0zba/NcYc8UY8xuwGajlnBKVVzEGxo+Hrl05U6IE5c+eJalCBXbs2KHh4AKOBMQOoKKIlBcRP6AHsNyqzVdAQxHJJyIFSLsEOeDcUlWed+NG2k1zR4xgx4MPUv7XX6nXvj1bt24lNDT0zu9XTnfHgDDGJANDgGjSvvRfGmN+EpFBIjIovc0B4L/Aj0AsEGmM2ZdzZas8Jz4eWreGzz5jTrlyhB85wqtvvsmSJUsoVKiQq6vzWmKMa7oCwsLCzM6dO13y2crNHD0KbdpgjhxhRNGiTIqPZ8aMGfTp08fVleVJIrLLGOPQ5BGdSalca9s26NCBpOvX6Zg/PztSU4mJiSEiIsLVlSl0sZZypagoaNqUiyLUvHKFE+XKERsbq+HgRvQMQuU46+3jh7WoRMdVn8PIkRwuVYrw06ep16YN8+bNo3Dhwq4uV2WgAaFylPX28ed/S4Bnn4W9G1hfsiRtTp9m6GuvMX78eCwWS9YHU7lOA0LlqIzbxwcmJvDvpR8QfnIfowsU4b1z55gWGUn//v1dXKWyRwNC5aibW7uV/z2OmYvepeSl8/Ty9Wd+aioxGzbQuHFjF1eosqIBoXJUqcAASv8Yy7SlH5CSkkRTk8qO+4pTq98HGg4eQANC5ahPkvZSY8HfOObnT+uka5wpV5vyXd5m5NN1XV2acoAGhMoZqakwahSPjhlDbJEgWl6KJ+XRdtTuNJQ3nqymOzx5CA0I5XzXrqWNVCxYwOLAQHolXOKjKVN48cUXXV2ZuksaEMqpVm/4gdDnnqbGif2M8PVncnIqK6OjadasmatLU9mgAaGcZv3ijdR4vifFE87TxceHrwoVJbTH37lcrKqrS1PZpFOtlXNs2ED4M23xuxJP49QUVoVWp0TvD0ktXPKOt7lT7ksDQt27yEhMq1acSE2lTtJ1Djzcivu7jsYScB/g2G3ulHvSgFDZl5oKw4fDCy+wLSCAejeu8ccTL1C0xWDE8ufVa2ABXxcWqe6FBoTKnqtXoWtX+Mc/mBUQQDvg//71BcXCOyKSeRvTy9eSWbY7zjV1qnuiAaHu3pkz0LgxZulShlksvFeiBFu3b2f0kGco6Hd7v3dSqtF+CA+loxjq7uzdi2nblqQzZ+hiDJciIvhu8WKKFy8OwKXEJJtv034Iz6RnEMpxa9ZgIiKIP3eOuklJBD/3HOvWrbsVDpC29sIWe88r96YBoYC0fRsixsVQfsQqIsbF3N5nMGUKpm1bDiUnU/PaNXr9859ERkbi5+eXqdmwlpUJ8M28r0OAr4VhLSvn9P+CygF6iaFu29Ql7mIiby7ZC0DHmiXg9ddh0iTW5c9PHx8fZixfTrt27Wwe6+Yai0w7SLWsrGsvPJTuaq2IGBdDnI0+ggoBsH7Xv2HFCiZbLHxYqhRfrVxJzZo1XVClchbd1VrdFVsdiCUSfuPjz0eTev4YQ4AfwsP5bulS7r9fb9zuTbQPQt3WgfjQ2SMs++JVyvx2nDbGcLl3b2JiYjQcvJAGhMrUsdjsl+9Y+J9hJF1LoH5qCo0++IDZs2eTP39+F1epXEEvMVRaB6IxHB85lqGrpvK9xUKXfPn4aOGXdOrUydXlKRfSgFCQnEzHyA9g1VS+slh4/YH7WbpyJbVr13Z1ZcrFNCC8XUICpls3JDqa8cDSRx5h81dfUbJkSVdXptyABoSXyXiXq9omgc8Xv0vBo78wCLjcvTtff/45AQE661Gl0U5KL3JzQlTcxURqnj7ItMkvYo4epqUxhL77LlFRURoOKhM9g/AiN+9y1frnrUxc+SFnUlNp52PhSqcRbBg1ytXlKTekZxBe5HT8VQZtX8TUr8axJzWZiID7+L3XeHiwvu7XoGzSgPAWSUl8HDOFEZtmMR9oFVwOS99J5C9ZCUD3a1A26SWGN4iPJ7VzZ9rv3Mh7wPiK9Sja9jV8/PxvNdH9GpQtGhB53NoV26jYrxshF07zPLAqojvFInohkvnkUfeNVLZoQORhm2d9xaMvPgPXr9Dcx4f9bV7hvmqP22zrokW9ys1pQORVUVGEP9+H4ykptAsoRMJToygYUsVuc3tbxSnvpgHh5jJObMq4+Yq95zEG8957yKhRfAd0K16GfF3fJX/hrFdi6pZwyhaHAkJEWgEfAxYg0hgzzk67x4DtQHdjzCKnVeml7O30tPP47yzeFXfb8z43rtP2XyPxmTuXL4A3KtfF78nX8PHL/OUXIOMVhW4Jp+y54zCniFiAT4HWQDWgp4hUs9NuPBDt7CK91c2JTRklJqUQ9d3J257PnxBP6W7t8Zk7l3eAn4YNY8rcLylYsFCmdgG+FnrVLUNIYAAChAQGMLZzDd0STtnkyBlEHeCwMeYogIjMBzoA+63aDQUWA485tUIvZm/oMcWqR7H873HMXPA3Siaco7fFwhORkTz77LMA+Pj46P6QKtscCYgQ4GSGx6eA8IwNRCQE6AQ0RQPCaUoFBtjcK9Iiciskwk/sZdri0STfSKSFf0HeX7uGhg0b3mrbsXaIBoLKNkdmUoqN56wHxSYBw40xKTba/nkgkQEislNEdp4/f97RGr2WvS3ke4aHEuBr4am9G5gz/23+dyORhkVL0WvhukzhoNS9cuQM4hQQmuFxaeC0VZswYH76PRmLA0+KSLIxZlnGRsaY6cB0SNvVOrtFewu7W8jXKkmv5dOpunoy64H+FR5l7Mz/8HRD+8OYSmWHIwGxA6goIuWBOKAH8HTGBsaY8jd/FpFZwErrcFDZc9slwrVr3OjWjaqLFxMJHBg6lKMffYTFYrF7DKWy644BYYxJFpEhpI1OWICZxpifRGRQ+uvTcrhGr2J3fgPAuXNca9UK/927GeHjw4NTp/LhgAGuLVjlaQ7NgzDGrAZWWz1nMxiMMc/ee1neKcs7XOW/RGKzZnDmDH0LFuTZFSto0qSJK8tVXkCXe7sRe/Mevv50HtfDwrh05gx9y5blb3v2aDioXKEB4UZszXvovmcN/5w5nEOJiTQvW50OX6ymQoUKLqhOeSMNCDeScT2EmFSGr5/O+OhPWW8MrWq1IKHb+4xZd0J3f1K5RgPCDSzbHZfpBrr+Sdf4dOHfeXHXcqYAfZ8YiG+rlxEfC4lJKbr7k8o1uprTxaw7JoMvxzN9/lvUunCSVyx+zOv8DgX+8kim9+juTyq3aEC4WMaOycrnfyUy6k2KJv5Bl4JBfNfzAwKKhd72Hl2arXKLBoSL3TwbaHRkB5OXvs8fKck8UaICcd3ewzfgvtva69JslZs0IFysVGAAj6+LYnRMJD8CPao25lqbV7BY0n41IhAY4MvFq0m6GlPlOg0IV0pJYdp3M6gRM5/lwOBGz+JT9ynS17SkMbB7ZAuXlai8mwaEq1y+TPyTT1JjyxYm5/NlbMe38H3w9pXy2t+gXEkDwhVOnSK+QQMKHz/OqGLF6LZpE6VvBGYazQDtb1CupwGRTVkuqspC6q5dXG7ShHx//MHwatUYvnEjwcHBPJT+uu7+pNyJBkQ2ZLmoKosv9PVFizA9enAxJYUZ7dvz/pdfkj9//luv6+5Pyt3oTMpssLeoyu4MR2O49Pe/49u1K3tTUljx9tuMXrYsUzgo5Y70DCIb7M1ktPl8cjLnn36a4IULWW6xYImKYnDXrjlcoVLOoQGRDfY2k71txCEhgf81bkyJPXuYet991N+0iVq1a+dSlUrdO73EyAZ7m8lmHHEwx49ztlIliu3Zw9jy5el06JCGg/I4egaRDbY2k21SJZgJ0Qd5ZcEemlw8zITPh5P/xnX+2bQpr6xahb+/v4urVuruaUBkU8YRh4yjGi1+WMuk6H9x1him9niO8fMiM8+MVMqDaEA4wYTogyTeSOb5mBm8s3M53yL0b/1Xgmq113BQHk0DwgnOXfiD9xb9nd6/7mZBPj9GdB+DKV2NRN23QXk4DYh7ZOLj+fzzl2hw4TRjCwYxpc9ELIWDAV1HoTyfBsQ9uH7gABfq1aPOpUsMKlWJ1T0+wOKb1hmp6yhUXqDDnNkUv3o1V2vVwv/SJeb27UvLFRsIDQ5CgJDAAMZ2rqHTppXH0zOIbDg5YQL3v/EGJ0U4OHEi/V55BYBOj5R2cWVKOZcGxN0whoN9+1J5zhy2+/rit2oVbZo3d3VVSuUYDQgbbC3l7lC1GAcaNqTazp2sDAqiVmwsoXoDG5XHaR+ElZuTnuIuJmJIW8o9bvYm9oaWpdrOnURVrUqTEyc0HJRX0ICwYr2UO/TMIeZ82o/Kv51jUqPm/OuZT6g+ZhMR42L0Dlcqz9NLDCsZl2w/um8D01dPItUYOtXvwbFGfUlMuA44vkmMUp5MzyCs3Jzc1PbrmUSt+ohzCO06DOdgo953t0mMUnmAnkFksGx3HInXbjDky5G8fux7Ynz9GfrMBFJDKpBiFQ436W3wVF7mlQFha5QCYOT8WEbOGEKX+NPMKhzM6L4fU7hYMd5t/xATog86tkmMUnmI1wWEvQ1ngxPOMP3fQ6l7/QqjSj/ErB7vI5Z8FMyf71Yfg25Lr7yN1/VB2Npw9v6fv2PW5Beodf0Kz9Z+ktm9xiPpt76Lu5hIxLgYAMZ2rkFIYIBOp1Zew+vOIKz7DMK+W0LkxplcQ+jcYjAHare+7T03zzLGdq7BNyOa5lapSrmc15xBLNsdR8S4GEyG5zqs/JCojTOJs+Sja+9x/Fqnrd3364iF8kZeERAZZ0cCkJLEX2f/Hx//9DVfBxSm24AZ/C+0Ok89GkJIFp2OOmKhvI1XBETGfge/y/FMmtKPV/53mMjiZRn00iyuFA4mKcXw9c/n+WZEU7shoSMWyts4FBAi0kpEDorIYREZYeP1XiLyY/p/20SklvNLzb6bf/mD4n5m7rT+tL96kbcqN2DMc5NJzed3WztHtrVXyhvcsZNSRCzAp0Bz4BSwQ0SWG2P2Z2h2DGhsjIkXkdbAdCA8Jwq25U430i0VGIDvlhV8seZjihnDsw2fYXP9Hrcd5+YZgq1t7fVGusobOTKKUQc4bIw5CiAi84EOwK2AMMZsy9B+O5BrO6c4ciPdZ/bNp9fqGfwhPnTu9BYHK9W77TjWZwh6I12lHLvECAFOZnh8Kv05e/oDa+6lqLuR1Y10U1JSWNqiBQPmzuB0gYIM+ut0DlWqR0hgAM/ULaNzGpS6A0fOIGzd2MHYeA4RaUJaQDSw8/oAYABAmTJlHCwxa3ZvpPu/86yqVIlOR4+yt2xZquzezYqgIKd8plLewpEziFNAaIbHpYHT1o1EpCYQCXQwxlywdSBjzHRjTJgxJiw4ODg79d7G1siC34UTTJ72HO2PHmXP449T48gRfDUclLprjgTEDqCiiJQXET+gB7A8YwMRKQMsAXobYw45v0z7rEccgn7ZzoKZQ2l1/Sr7X3qJh7/+GiyWLI6glLLnjpcYxphkERkCRAMWYKYx5icRGZT++jRgJFAMmJJ+q7lkY0xYzpWdeeSiSIAv/r4+BK2fw5wt/6GICHHTplFt4MCcLEGpPE+MsdmdkOPCwsLMzp07s/Ve65ELY1JpsGIC0w9s4XL+/ORfu5bARo2cWa5SeYaI7HL0D7hHLtbKOHKRej2RnnOHMfb8r/xQoAgP7f8Bv7JlXVyhUnmDR061vjlyYS6e5e1/92f8+V9ZFVyOngMjNRyUciKPDIhSgQH4Hv+R6ZGDGJiYwORK9Xm53ycUu7+oq0tTKk/xyEuMJy9+S/v5f6MKhtcjerKoQS9dK6FUDvCogDDGEDlwIP1nzKCAj4X/6/E2q0PrEKJrJZTKER4TEImJiUxt1oxB27ZxpVAh/Lds4dOHH3Z1WUrlaR4RECY1ldkPPcSrx44RV6YMpWJjkQcecHVZSuV5HhEQIkKr0FBOhYRQeu1aCNCNW5TKDR4xirFsz2n6tHibRvXfIOLjb/WemErlErefSWk9axLA10co5J+Pi1eTdDMXpe5SnppJaWu/h6RUQ/zVJEBvoqtUTnL7SwxHdpLWLemVyhluHxCO7iStW9Ir5XxuHxC2dpi2RbekV8r53L4PwnqH6SIBvly5kUxSyp+dqzrNWqmc4fYBAbfvMH2nbe6VUs7hEQFhTbekVyp3uH0fhFLKdTQglFJ2aUAopezSgFBK2aUBoZSySwNCKWWXBoRSyi4NCKWUXRoQSim7NCCUUnZpQCil7PLItRhKOYOrF/25+vMdoQGhvJL1XqcEYbrmAAAEcUlEQVS5vXWhqz/fUXqJobySrb1Oc3PrQld/vqM0IJRXsrdFYW5tXejqz3eUBoTySva2KMytrQtd/fmO0oBQXsnWXqe5uXWhqz/fUdpJqbyS9V6nuT2K4OrPd5Tb31lLKeVcd3NnLb3EUErZpQGhlLJLA0IpZZdDASEirUTkoIgcFpERNl4XEfkk/fUfReQR55eqlMptdwwIEbEAnwKtgWpATxGpZtWsNVAx/b8BwFQn16mUcgFHziDqAIeNMUeNMTeA+UAHqzYdgC9Mmu1AoIiUdHKtSqlc5khAhAAnMzw+lf7c3bZRSnkYRwJCbDxnPXnCkTaIyAAR2SkiO8+fP+9IfUopF3JkJuUpIDTD49LA6Wy0wRgzHZgOICLnReT4XVULxYHf7vI9uUnruzda371xtL6yjh7QkYDYAVQUkfJAHNADeNqqzXJgiIjMB8KBS8aYM1kd1BgT7GiRN4nITkdngLmC1ndvtL57kxP13TEgjDHJIjIEiAYswExjzE8iMij99WnAauBJ4DBwFejnzCKVUq7h0GItY8xq0kIg43PTMvxsgMHOLU0p5WqeNpNyuqsLuAOt795offfG6fW5bDWnUsr9edoZhFIqF7llQLj72g8H6uuVXtePIrJNRGq5U30Z2j0mIiki0sXd6hORx0Vkj4j8JCKb3KU2ESkiIitE5If02nK1Q15EZorIORHZZ+d15343jDFu9R9pIyVHgL8AfsAPQDWrNk8Ca0iboFUX+M7N6qsPBKX/3Nrd6svQLoa0zucu7lQfEAjsB8qkP77fjWp7Cxif/nMw8Dvgl4v/fo2AR4B9dl536nfDHc8g3H3txx3rM8ZsM8bEpz/cTtrEsdziyL8fwFBgMXAuF2sDx+p7GlhijDkBYIzJrRodqc0A94mIAIVIC4jkXKoPY8zm9M+0x6nfDXcMCHdf+3G3n92ftETPLXesT0RCgE7ANHKfI/9+lYAgEdkoIrtEpI8b1TYZqEraTOG9wF+NMam5U55DnPrdcMdNa5229iOHOPzZItKEtIBokKMVWX2sjees65sEDDfGpKT9IcxVjtSXD3gUeAIIAL4Vke3GmENuUFtLYA/QFHgQWCciW4wxCTlcm6Oc+t1wx4Bw2tqPHOLQZ4tITSASaG2MuZBLtYFj9YUB89PDoTjwpIgkG2OWuUl9p4DfjDFXgCsishmoBeR0QDhSWz9gnEm74D8sIseAKkBsDtfmKOd+N3Krc+UuOmHyAUeB8vzZUfSQVZs2ZO6IiXWz+sqQNu28vjv++1m1n0XudlI68u9XFdiQ3rYAsA+o7ia1TQXeTf/5AdLWJxXP5d9xOex3Ujr1u+F2ZxDGzdd+OFjfSKAYMCX9r3SyyaVFPg7W5zKO1GeMOSAi/wV+BFKBSGOMzWG93K4NeA+YJSJ7SfsSDjfG5NoKTxGJAh4HiovIKWAU4JuhPqd+N3QmpVLKLnccxVBKuQkNCKWUXRoQSim7NCCUUnZpQCil7NKAUErZpQGhlLJLA0IpZdf/A2LUUGV77LhzAAAAAElFTkSuQmCC\n",
            "text/plain": [
              "<Figure size 288x288 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    }
  ]
}